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Mach number dependence of the onset of dynamo action 
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ABSTRACT 

The effect of compressibility on the onset of nonhelical turbulent dynamo action is 
investigated using both direct simulations as well as simulations with shock-capturing 
viscosities, keeping however the regular magnetic diffusivity. It is found that the criti- 
cal magnetic Reynolds number increases from about 35 in the subsonic regime to about 
70 in the supersonic regime. Although the shock structures are sharper in the high 
resolution direct simulations compared to the low resolution shock-capturing simula- 
tions, the magnetic field looks roughly similar in both cases and does not show shock 
structures. Similarly, the onset of dynamo action is not significantly affected by the 
shock-capturing viscosity. 



1 INTRODUCTION 



Transonic and supersonic turbulence is widespread in many 
astrophysical settings where magnetic fields are believed 
to be generated and maintained by dynamo action. Ex- 
amples include supernova- dr i ven turbulence in the inter- 
stellar medium jKo^^ t alJ Il999l : Ide Avillez &; Mac Lowl 
2002 :lBa lsara et al.ll2004T) and turbulence in galaxy clusters 
(^Det^TOr^t^dT^Q^ where the driving comes mainly from 
cluster mergers. In these cases the root mean square (rms) 
Mach number is of the order of unity, but in the cooler parts 
of the interstellar medium the sound speed is low and the 
flows can therefore easily become highly supersonic and may 
reach rms Mach numbers of around 20 (jPadoan et al.l ll997: 
IPadoan fc NordlundlEoO^l . 

Simulations do not give a clear picture of how the 
excitation conditions for dynamo action change in the 
highly supersonic regime. While there are cl ear examples 
of dy namo action in supersonic turbulence (Bal sara et alJ 
2004), the dynamo seems to be less efficient i n the more 
strongly supersonic regime (jPadoan et al-lEooih . It is plau- 
sible that newly generated magnetic field gets too quickly 
entrained by the shocks w here the field is then dissipated 
iPadonTi fc NorriliiTKiiflflft. 

Supersonic turbulence has a significant irrotational 
component. Purely irrotational turbulence is also referred 
to as acoustic turbulence a nd can be described us- 
ing weak t urbulence theory (Zak harov &; Sagdeevl Il97(l 
lBvkovl H988 > ). Since supersonic turbulence contains shocks, 
such flows can also be descr ibed as shock turbulence 
(|Kadomtsev &; Petviashvilil ll973F) . If these flows were de- 
scribed as predominantly irrotational, the growth rate 
would increase with th e Mach number to the fo urth power 
(iKazantsev et al J 1198.4 iMoss fc Shukurovlll99fft . However, 
this tendency has never been seen in simulations. 

There are however caveats to both the analytical and 
the numerical approaches. Firstly, it is clear that even highly 



supersonic turbulence is not fully irrotational and that still 
70-80% of the kinetic energy comes from the solenoida l 
component dPorter et alJll998l : IPadoan fc Nordiundlll999T) . 
Therefore purely acoustic turbulence cannot be used as an 
approximation to supersonic turbulence. The other problem 
that we shall be concerned with here is the question to what 
extent does the numerical shock viscosity used in many sim- 
ulations affect the conclusion regarding dynamo action. 

In order to address the latter question we perform direct 
simulations at sufficiently high resolution so that no numer- 
ical shock viscosity is needed. Here, the term 'direct' means 
that one uses the microscopic viscosity, i.e. one ignores the 
fact that in many astrophysical applications one will never 
be able to reach realistic Reynolds numbers. We also com- 
pare with calculations where a shock-capturing viscosity is 
included. Here the viscosity is locally enhanced in a shock, 
allowing less diffusion between the shocks. It turns out that 
both direct and shock-capturing viscosity simulations pre- 
dict an approximately similar increase in the critical mag- 
netic Reynolds number for dynamo action. In these cases the 
magnetic diffusivity is the same as in the direct simulations. 



2 METHOD 

We consider an isothermal gas with constant sound speed 
c s . The continuity equation is written in terms of the loga- 
rithmic density, 

Dlnp 



Bt 



-V • u, 



(i) 



and the induction equation is solved in terms of the magnetic 
vector potential A, where B — V x A, and 

OA 



dt 



= ux B + 77 V A, 



(2) 



77 being the magnetic diffusivity. The momentum equation 
is solved in the form 
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-c s 2 V \np + p 1 (J x B + F visc + /) : 



(3) 



Du 

where D/Dt = d/dt + u • V is the advective derivative, 
J = V x B//i the current density, I? the magnetic field, 
po the vacuum permeability, .Fvisc is the viscous force (see 
below), and / is a random forcing function with 



f(x, t) = Re{Nf k(t) exp[ife(t) • x + i0(t)]}, 



(4) 



where cc is the position vector. The wavevector k(t) and the 
random phase — tt < <fi(t) ^ tt change at every time step, so 
f(x,t) is, up to discretization errors, ^-correlated in time. 
We force the system with nonhelical transversal waves, 



f k = (kx e)/y/k 2 -(k-e) 



(5) 



where e is an arbitrary unit vector not aligned with k] note 
that \fk\ 2 = 1. On dimensional grounds the normalization 
factor N is chosen to be N = /oc s (|fc|c s /(5£) 1//2 , where fo 
is a nondimensional forcing amplitude. We use the Pencil 
Code, 1 which is a high-order finite- difference code (sixth 
order in space and third order in time) for solving the com- 
pressible hydromagnetic equations. 

In the direct simulations with constant viscosity, zy, the 
viscous force per unit mass is given by 



^visc = pv iy 2u + 1 vv • u + 2S • V M 



where p is the density and 
du 



duj 



+ 



dxi 



(6) 



is the traceless rate of strain matrix. In runs with shock- 
capturing viscosity we simply add to this a spatially 
dependent diffusion term such that the effective viscos- 
ity is enhanced only in th e neighbourhood of a shock 
( Ric htmver &; Morton! Il967lh This technique artificially 
broadens the shocks such that they can be resolved numeri- 
cally, and hence all the conservation laws are obeyed, which 
is important for satisfying the right jump conditions. It is 
sufficient to broaden the shocks using only a bulk viscosity, 
£, rather than a locally enhanced shear viscosity. Thus we 
write 



-, (shock) 



F v l £+pCVV-u + (V-u)V(pC), 
so the full stress tensor is given by 
nj — 2pvSij + pCSijV • u. 



(8) 



(9) 



Here C is the sho ck viscosity. Following 
iNordlund &; Galsgaardl (|l995l) . we assume ( to be pro- 
portional to the smoothed (over 3 zones) maximum (over 5 
zones) of the positive part of the negative divergence of the 
velocity, i.e. 



C = c s hock ^max[(-V • u) + ] 



(10) 



where c s hock is a constant defining the strength of the 

shock viscosity. Th is is also the technique used by 

IPadoan fc Nordlundl fcOOfl . 

The simulations are governed by three important pa- 
rameters, the Mach number Ma = u rms /c s , the Reynolds 
number Re = u rm s/ {yki ), and the magnetic Prandtl number 

1 http : //www . nordita . dk/data/brandenb/pencil- code 



Table 1. The Mach and Reynolds numbers, as well as the growth 
rate (in units of Ao = itrms&f) for runs with different forcing 
strength (fo), viscosity (y) and c s hock, an d number of meshpoints 
N . All runs have Ptm — 1- The runs with only direct viscosity 
are Runs la-3a. 



Run 


N 


fo 


V 


Ma 


Re 


A/Ao 


^shock 


la 


512 


0.2 


0.006 


0.72 


78 


+0.030 




lb 


128 


0.2 


0.006 


0.70 


78 


+0.010 


0.8 


lc 


128 


0.2 


0.006 


0.68 


74 


+0.014 


3 


2a 


256 


0.2 


0.01 


0.65 


43 


+0.006 




2b 


64 


0.2 


0.01 


0.62 


41 


-0.004 


0.8 


2c 


64 


0.2 


0.01 


0.65 


43 


-0.005 


3 


3a 


512 


0.5 


0.01 


1.11 


75 


+0.005 




3b 


128 


0.5 


0.01 


1.14 


76 


+0.010 


0.8 


3c 


64 


0.5 


0.01 


1.10 


74 


+0.002 


3 



PrM = v/rj. The magnetic Reynolds number is defined as 
ReM = Re PrM, and the critical magnetic Reynolds number 
is the value above which a weak seed magnetic field grows 
exponentially. 

The Mach number can be increased either by increas- 
ing the strength of the forcing or by lowering the sound 
speed. We choose the former and keep the sound speed con- 
stant. Moreover, we use c s as our velocity unit. The lowest 
wavenumber in the box, ki = 2tv/L, is used as our inverse 
length unit. This implies that time is measured in units of 
(c s /ci) _1 . Density is measured in units of the mean density, 
which is also equal to the initial density po. The extent of 
the box is in each direction L — 2tt, so the smallest possi- 
ble wavenumber in the box is therefore k — 1. The maxi- 
mum possible wavenumber depends on the resolution and is 
k = 256 for our largest run with 512 3 mes hpoints. Larger res - 
olution is technically already possible ([Haugen et al.l l2003 N ). 
but still too demanding in terms of computing time if one 
wants to cover many turnover times. In all our simulations 
the flow is forced in a band of wavenumbers between 1 and 2. 
As initial conditions we have used zero velocity and a weak 



random magnetic field with (-B 2 ) /(popoc^) 
periodic boundaries in all three directions. 



10"*. We use 



3 RESULTS 

In TableQwe give the Mach and Reynolds numbers as well as 
the growth rate for runs with three different combinations 
of fo and is, each with and without added shock viscos- 
ity. The run parameters are generally chosen such that the 
growth rates are close to zero so that the critical magnetic 
Reynolds number for dynamo action can accurately be de- 
termined via interpolation. It turns out that the addition of 
shock viscosity reduces the normalized growth rate of the dy- 
namo, A/Ao, only by about 0.003-0.016 (see Tabled, which 
is about 3-16% of the typical normalized growth rates of 
abou t 0.1 in the more s trongly supercritical cases (see Fig. 3 
of lHaugen et all |2004 ) . Here, Ao = w r msb) is the typical 
stretching rate. The reduction of the growth rates can, at 
least partially, be explained by a reduction in ReM due to 
reduced rms velocity when shock viscosity is added. 

We make a more detailed comparison of two equivalent 
simulations, one with and one without shock viscosity. In 
each of the two cases the Mach number is the same, e.g. 
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Figure 2. Value of V • u along a line through the box whose 
position is indicated in Fig.Q] Here Ma = 1.1 and we are compar- 
ing constant viscosity (solid line) and shock-capturing viscosity 
(dashed line). 
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Figure 3. Same as Fig. [I] but for B z . Left figure is from a run 
with 512 3 meshpoints and normal viscosity while right figure is 
from a run with 128 3 meshpoints and shock viscosity. Since this is 
in the linear regime we have normalized by the rms value. There 
is surprisingly little difference between the two. 



Figure 1. Grey scale (colour scale in electronic version) rep- 
resentation of V • u in an xy cross-section through z = for 
Ma = 1.1 using constant viscosity (upper panel: Run 3a, 512 3 
meshpoints) and shock-capturing viscosity (lower panel: Run 3b, 
128 3 meshpoints). Both simulations have been run from an initial 
state with zero velocity, to a state where the kinetic energy has 
saturated, but the magnetic energy is still in the linear regime 
(i.e. there is no back reaction on the velocity). Note the drastic 
difference between the two regarding the sharpness of regions of 
strong convergence. In both panels, the dotted line indicates the 
position of the scan shown in Fig. |2] 



Ma ~ 1.1, fo = 0.5. Corresponding cross-sections of V • it 
are shown in Fig. Q In the direct simulations the regions 
with V • u < are quite sharp, compared with the shock- 
capturing simulations where they are smoother. This is seen 
more clearly in Fig. |^1 where we plot V • u along a cross- 
section (shown as a dashed line in Fig.QJ. Since the runs are 
turbulent, however, direct comparison between individual 
structures is meaningless and one can only make qualitative 
comparisons. 

The number of strong convergence regions and shocks 
is roughly similar in both cases, but this is only because we 
have used the same grey/colour scale in both plots. In the 
direct simulation, the dynamical range is —62 ^ V u ^ 4.2, 



which exceeds the range of the grey/colour scale, whereas in 
the shock- capturing simulation the dynamical range is only 
— 17 ^ V u ^ 4.2. If one were to compare the two cases 
such that in each the grey/colour scale is exactly within 
the dynamical range of that simulation, one would see more 
structures in the shock-capturing simulation. 

The dynamo-generated magnetic fields are generally 
rather smooth in both cases and in that respect rather simi- 
lar to each other; see Fig.|3] Indeed, the filling factors where 
the field exceeds its rms value are similar in the two cases 
(0.09 in the left hand panel for the direct simulation simula- 
tion and 0.13 in the right hand panel for the shock capturing 
simulation). The probability density functions of the three 
components of the magnetic field are in both cases stretched 
expo nentials (Fig. El , which is in agreement with earlier re- 
sults (Bra ndenburg et aTlll996l) . 

In order to elucidate the reason for the almost com- 
plete absence of any shock-like structures in the dynamo- 
generated magnetic field, we show in Fig. |5] the different 
contributions to the right hand side of dB z /dt for Run 3a 
with normal viscosity. The grey /colour scale is the same in 
all panels. It is evident that the dominant term in the induc- 
tion equation is the advection term, —u -V5 Z . The stretch- 
ing term, B • Vii 2 , gives the weakest contribution. The com- 
pression term, — S z V-w, gives an intermediate contribution, 
but, more importantly, there are only a few barely notice- 
able contributions from shocks. This is mainly because in 
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Figure 4. Probability distribution function of B x , B y and B z , 
normalized by £> rms , for Run 3a (upper panel) and Run 3b (lower 
panel). 
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Figure 5. The different terms in the induction equation for 
Run 3a (with 512 3 mesh points and constant viscosity). We clearly 
see that the dominant term is the advection term. 




Figure 6. Same as Fig.0 but for J z . Since this is in the linear 
regime we have normalized by the rms value. The current sheets 
are similar in thickness, but there are more of them in the shock- 
capturing simulation. 



the locations of the shocks (for example at x ~ — 1.5 ?/ ~ 2; 
see upper panel of Fig. the normal component of the 
magnetic field is weak (Fig. 0, making the shocks less pro- 
nounced in the product of the two (i.e. in i^V -u). The end 
result is a relatively weak contribution to the right hand side 
of dB z /dt. 

In the cross- sect ions of the current density we see very 
similar structures in the simulations with normal and shock 
capturing viscosities; see Fig. This is supported by the 
filling factors which are ^0.1 for both the left and the right 
hand sides of the figure. 

In Fig. [7| we plot the critical magnetic Reynolds num- 
ber ReM,crit as a function of Mach number. It turns out 
that the critical magnetic Reynolds number increases from 
about 35 in the subsonic regime to about 70 in the super- 
sonic regime. Whether or not the critical magnetic Reynolds 
number increases even further for larger Mach numbers is 
unclear, because larger resolution is required to settle this 
question for PrM = 1. However, in the direct simulations 
with PrM = 5 (Fig. 0, and in the shock-capturing simula- 
tions with PrM = 1 the critical magnetic Reynolds number 
is roughly unchanged when Ma is increased beyond Ma ~ 1. 
We refer to this possibility of having two distinct values of 
the critical magnetic Reynolds number for subsonic and su- 
personic turbulence, with a reasonably sharp transition at 
Ma ~ 1, as 'bimodal' behaviour. 

Comparing direct and shock-capturing simulations we 
find the general appearance of the cross-sections of B re- 
markably similar. The simulations with shock-capturing 
viscosity have, however, slightly larger critical magnetic 
Reynolds numbers. This is probably explained by the smaller 
velocity gradients and hence smaller stretching rates in the 
simulations with additional shock-capturing viscosity. Nev- 
ertheless, both direct and shock-capturing simulations show 
a roughly similar functional dependence of the critical mag- 
netic Reynolds number on the Mach number. This suggests 
that shock-capturing simulations provide a reasonable ap- 
proximation to the much more expensive direct simulations- 
at least as far as the onset of turbulent dynamo action is 
concerned. 

Once the Mach number becomes comparable to unity, 
one needs very high resolution to resolve the rather thin 
shock structures. However, as seen above, the magnetic field 
structures are by far not as thin as the hydrodynamic shocks. 
This implies that one can decrease 77 to values far below 
v before the magnetic field and current structures would 
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Figure 7. Critical magnetic Reynolds number ReM,crit as a func- 
tion of Ma for simulations with Pr^f = 1. Note that Rejucrit 

de- 
pends strongly on Mach number for Ma ~ 1. The simulations with 
shock-capturing viscosity give approximately the correct growth 
rates. The simulations that provide these data points have reso- 
lutions ranging from 64 3 to 512 3 mesh points. (Some of them are 
listed in Table □) 
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Figure 9. Energy spectra for Runs la and lc. The dotted lines 
give the result using shock-capturing viscosity (numbers in paren- 
theses in the table). 
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Figure 8. Same as Fig. |7|but for simulations with Ptm = 5. 

become unresolved. In other words, with given resolution, 
one can reach a magnetic Reynolds number that is much 
larger than the kinematic Reynolds number. Alternatively, 
if one only wants to reach just weakly supercritical values of 
ReM for dynamo action (which turn out to be roughly below 
one hundred), Re can be much smaller and just a few tens. 
This regime corresponds to large values of Ptm- For Ptm = 
5, for example, we are able to reach Mach numbers of up 
to 2.1. The resulting plot of ReM,crit versus Ma is shown in 
Fig. El This plot confirms again the approximately 'bimodal' 
behaviour of the critical magnetic Reynolds number, but 
now ReM,crit ~ 25 in the subsonic regime and ReM,crit ~ 50 
in the supersonic regime. 

The kinetic and magnetic energy spectra show close 
agreement between the constant and shock- capturing vis- 
cosity solutions at low wavenumbers; see Fig. E] The trun- 
cation at higher wavenumbers is presumably a direct con- 
sequence of the locally enhanced shock viscosity which at- 
tempts to increase the width of the shocks. Furthermore, 
since shocks hardly manifest themselves in the magnetic 
field, the agreement between the magnetic energy spectra 
in the direct and shock-capturing simulations extends to 



even larger wavenumbers. As discussed in the beginning of 
Sect. El the direct simulations require vastly more resolu- 
tion, but once the resolution is sufficient, either of the two 
spectra agrees well when dou bling the reso l ution ; see also 
the middle panel of Fig. 13 in lHaugen et al.1 {2004). 



4 DISCUSSION 

The present results have shown that in both direct and 
shock-capturing simulations the onset of dynamo action re- 
quires a somewhat larger magnetic Reynolds number when 
the Mach number exceeds a critical value around unity. In 
the subsonic regime the critical magnetic Reynolds number 
for nonhelical dynamo action is around 35, but it increases 
to about 70 in the supersonic regime. Once the Mach num- 
ber exceeds unity, the critical magnetic Reynolds number no 
longer seems to depend on the Mach number. This confirms 
the notion that in the supersonic regime dynamos experience 
an additional sink. This additional sink is plausibly related 
to the sweeping up of magnetic field by shocks, but this ef- 
fect does not gain in its importance with increasing Mach 
number once we are already in the supersonic regime. 

The relative importance of the adverse effects of shocks 
can be associated with a Reynolds number-like quantity 

R divu = ((V-u) 2 ) 1/2 /(^ f 2 ), (11) 

whose dependence on the actual Reynolds number is shown 
m Fig. HJ3 Here we also compare with a similar number 
quantifying the relative importance of vortical motions, 

R» = (u 2 ) 1/2 / (./*?) . (12) 

Both numbers seem to approach a limiting behaviour pro- 
portional to Re 1-2 . On theoretical grounds, one would have 
expected a slope of 3/2, because the ratio of rms velocity to 
rms vorticity is proportional to the Taylor microscale which, 
in turn, is known to diminish proportional to Re -1 / 2 . The 
departure from this expectation might indicate that we are 
not yet in the asymptotic regime with a well developed in- 
ertial range. This is also plausible from Fig. showing that 
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Figure 10. Dependence of i^div-u and Ruj on Re. Note that both 
quantities remain well separated as they tend to approach a slope 
of 1.2. 



most of the spectrum is dominated by a rather extended 
dissipative subrange with only a rather short inertial range. 
Thus, it might not be surprising that the Taylor microscale 
does not yet show a clear asymptotic Re -1 / 2 scaling. 

Figure fTTTl suggests that RdWu remains smaller than R u 
by a constant factor of around 1.5. Thus, ((V • it) 2 ) / (u; 2 ) ~ 
0.44. This ratio would be exactly 1/2 if the mean square val- 
ues of longitudinal and transversal velocity derivatives were 
equal, i.e. (u 2 ,^) = (u 2 ^). Here we have assumed isotropy 
and that mixed terms cancel, which implies ((V • u) 2 ) ~ 
^( u x,x) an d (u? 2 ) ~ 6(^ 2 >y ), giving a ratio of 1/2. 

Obviously, it will be important to confirm the asymp- 
totic behaviour of R 2 J /R 2 lWu at larger Reynolds and Mach 
numbers and to relate this to the value of the critical 
magnetic Reynolds number for dynamo action. At the mo- 
ment our simulations are simply limited by the resolution 
(512 3 meshpoints) which in turn is limited by the com- 
puting power available on modestly big supercomputers. 
Even though somewhat larger resolutions are already pos- 
sible (1024 3 meshpoints), those runs remain prohibitively 
expensive if one needs to run for many turnover times. 

Although we have here focused on the critical value of 
the magnetic Reynolds number for dynamo action, rather 
than the growth rate, we can see that there is currently no 
evidence for an increasing growth rate with increasing Mach 
number. From Runs la and 3a in Tabled one sees that, as 
Ma increases from 0.72 to 1.1, the growth rate actually de- 
creases, even though the magnetic Reynolds number is the 
same. By comparison, early wor k on small scale dynam o ac- 
tion in irrotational turbulence (Ka zantsev et all 11985) has 
shown that the growth rate should increase with the Mach 
number to the fourth power. Although we cannot confirm 
this, it is important to bear in mind that irrotational turbu- 
lence is not likely to correspond to the limit of large Mach 
numbers. Instead, as discussed above, the irrotational contri- 
bution to the flow velocity may converge to a finite fraction 
of the total flow velocity in the limit of large Mach numbers. 



5 CONCLUSIONS 

The present simulations suggest that the critical magnetic 
Reynolds number for the onset of dynamo action switches 
from about 35 in the subsonic regime to about 70 in the 
supersonic regime, and stays then approximately constant 



(at least in the range 1.0 ^ Ma ^ 2.6 that we were able 
to simulate). It turns out that the simulations with shock- 
capturing viscosity yield almost the same critical magnetic 
Reynolds number as the direct simulations. This is partly 
explained by the fact that the shock structures are approxi- 
mately equally smooth for both direct and shock-capturing 
viscosities. Finally, we note that there is as yet no support 
for the expectation that the growth rate is proportional to 
the Mach number to the fourth power, as was suggested pre- 
viously bv_jiialv_tic__th^rv^ 

tional (|Kazantsev et ajjEjj |Moss & Shukurov|[l99i) . It is 
suggested that this is because supersonic turbulence cannot 
accurately be described as being fully irrotational. Indeed, 
we find that the Reynolds number based on the vortical 
component of the flow is always larger than the Reynolds 
number based on the irrotational component. 

In order to verify the proposed asymptotic indepen- 
dence of the critical magnetic Reynolds number on the 
Mach number, it would be useful to extend both the di- 
rect and the shock-capturing simulations to larger values of 
the Mach number. An obvious reason why this may not have 
been done yet is that when shock capturing viscosities are 
used, one normally uses at the same time also a similarly 
defined shock-capturing magnetic diffusivity. Consequently, 
a definition of the usual magnetic Reynolds number is no 
longer possible, and the connection between effective mag- 
netic Reynolds number (based on the shock-capturing mag- 
netic diffusivity) and the actual one is unclear. This is why 
in the present work we kept the magnetic diffusivity and 
the kinematic viscosity equal to the microscopic one. Even 
when we compared with shock-capturing simulations, the 
magnetic diffusivity was still kept constant. In conclusion, 
we feel that the use of shock-capturing viscosities in dynamo 
simulations with constant magnetic diffusivities provides a 
reasonable tool for investigating supersonic hydromagnetic 
turbulence. 
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